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Abstract. The evolution of the energy spectrum of cosmic rays accelerated 
by the first order Fermi mechanism, by a supernova remnant shock wave, 
including adiabatic deceleration effects behind the front, is carried out 
by means of a time-dependent numerical code. The calculations apply to 
the adiabatic stage (or Sedov stage) of the supernova explosion, and the 
energetic particle spectrum is calculated in the test particle limit 
(i.e., the back reaction of the cosmic rays on the flow is not included). 
The particles are injected mono-energetically at the shock. We show the 
radial distribution, and the spectrum of the accelerated and decelerated 
particles. 


I. Introduction. Diffusive acceleration of cosmic rays at collisionless 
astrophysical shock waves, where the particles pass (diffuse) across a 
plane shock repeatedly is a promising candidate for the origin of cosmic 
rays (Axford, Leer and Skadron, 1977; Krymskii, 1977; Bell, 1978; 
Blandford and Ostriker, 1978). ’ * 

Supernova explosions have long been supposed to be the origin of 
cosmic rays, primarily because of the energy they liberate. The 
diffusive shock acceleration mechanism has been applied to spherically 
symmetric supernova shocks (Krymskii and Petukhov, 1980, Prishchep and 
Ptuskin, 1981; Drury, 1983) and analytical solutions can be obtained with 
the assumption that K/RR is a constant (see also Blandford and Ostriker, 
1980 Bogdan and Vblk, 1983; Moraal and Axford, 1983). 

In this paper, we solve the time dependent problem of cosmic rays 
diffusively accelerated in a spherical shock by numerical integration. 


II. The Model. Under the assumption that there is sufficient scattering 
that the pitch-angle distribution is nearly isotropic (diffusion limit) 
the cosmic ray transport equation (Parker, 1965) may be written (for 
spherical symmetry) as : _ c . / „ £ \ . . _ . 

It = * S , u> 


where i 8 the phase space distribution function, U is the flow 

velocity, and S is the source. 

If there is a shock at radius R, we require f to be continuous at 
the shock (Toptygin, 1980), and the jump in streaming flux can be 
obtained by integrating (1) across the shock. Subscript 1 corresponds to 
the region outside the shock (upstream region), while subscript 2 
corresponds to the region inside the shock (downstream region), the 
boundary conditions at the shock are: 



III Supernova Shock If the shock front is moving, i.e., R = R(t), then 
it is useful to use the following dimensionless quantities: 


l-Vm) , X = and ^U>(t>/P.) , (3) 

where K p and P 0 are constantsto be defined later. This transformation 
breaks down if itR/R=l because the Jacobian of the transformation is then 
zero. Note that if K is constant and so is RR/K, then R oC, t* 1 ^ . 

Assume a velocity profile TJ= RV(\)8(1- O (e.g. Sedov's blast wave; 
Sedov, 1959) and consider particle injection at the shock, i.e., 
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S =QUt-K^Cp-FO = -^r fcCVDU*) , 


(4) 


where the units of Q are #/ sec/ length^/ momentum^. If VtOgoes to zero 
faster than > k >1 (which is true fair Sedov's blast wave 

approximation), then V(o)=:o *v>A o (i.e., the boundary condition 

at is ). If we assume self-similar solutions to the shock, R 


= At 


then by equation (3) we have: 


ty - 

L ~ A x 


, R , A (^)"‘ / "-“ > Vr=^^). <5) 

Substituting (3), (4) and (5), together with K-k x e x \ 'V WVl) 

into ( 1 ) and ( 2 ), the equations for the transport of cosmic rays in a 
moving spherical shock can be obtained. In order to have a forward in 
time equation, ©< has to be smaller than 1/2 (e.g. Sedov's blast wave, 
*<= 2/5). 

IV, Results and Discussion . We employ i finite differencing scheme 
implicit in the spatial variable but explicit in the momentum variable; 
it is also second order in spatial variable and first order in the 
momentum variable. To satisfy causality, we use downstream differencing 
for momentum except at the shock, where we use upstream differencing. We 
replace SC^) by a Gaussian with a spread of the order of the step 
size in q during computation. 

In the case of a moving spherical shock, we use the supernova shock 
as an example. We concentrate on the adiabatic stage (or Sedov's stage) 
of a supernova explosion only. 

The self-similar solution (Sedov, 1959) to the blast wave equations 
is R « At* where A - (E/f, )'* and <* = 2/5. The energy, E, is 
proportional to the explosion energy, E 0 . If T = 5/3, then E « 2.02 Eo 
. p, is the density of the undisturbed (i.e. upstream) gas. The 
velocity Vi%) is given by V(i) “ 5*/*. % v(%) where v(^) is a complicated but 
monotonically increasing function and £ V(A) ^ . 

The adiabatic stage of a supernova remnant begins roughly when the 
mass of interstellar gas swept by the shock is equal to the expelled mass 
of the supernova (e.g. Spitzer, 1978). At the beginning of this stage 
the shock radius is: w 

r. - / sgMoy /i ( 6 ) 

v- 4--TC P, ' 

where (5» is the expelled mass in solar mass units (Mq). Substituting (4) 
and r. “/U* into ( 6 ), we get: 

The adiabatic stage will last roughly until the temperature behind the 
shock drops below the recombination temperature (i.e. the shock becomes 
radiative) (e.g. Spitzer, 1978). Let subscript f represent this. The 
relation between upstream and downstream temperature for a strong shock 
is (e.g. Landau and Lifshitz, 1959): 

l YONO Pjj£ 

Yp, • (8; 


T*. _ iYCY-Q, vi- „ 

T, C Y+ ()*" J 

Assume an ideal gas law in the undisturbed (upstream) region 
where /a is mean molecular weight. This, together with R^.= <* 


and T£ = T r (the recombination temperature), gives: 


r - [< f ( J- QtJl \ 
T ( “ e/ V W feTr/ 


*4 


(9) 




Both % and of the adiabatic stage are proportional to P, * . We have 
studied cases with Y- 5 A (which gives E = 2.02 E 0 (Sedov, 1959)) and the 
ratio of the density of helium to that of hydrogen^f^/p* f = 0.4 . 
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la Figure 1 to 3 we used E 0 ■ 3 x 10^ erg, p B 0*3, Kj = K 2 = K a ” 
cm 2 sec" 1 , ai « a 2 s 0, Tr=10 2b K (the recombination temperature oT proton 
and electron; and Z the number density of hydrogen n„ = 3 x 10" J cm"- 3 (hot 
interstellar medium); i*e*, x loT* and r 5 .*>6xt(f 5 • An absorb- 

ing boundary is set up at a radius 2R. 

By studying Figure 1, we see that only a very small portion of 
particle (high or low energies) stays outside the shock, the absorbing 
boundary can be considered as being at infinity where the particle 
distribution, f, is supposed to be zero* This can be understood by the 
fact that there is no convection outside the shock and the diffusion time 
scale (R*7K P ) is much larger than the accelerating time scale (K 0 /U ); 
and the adiabatic stage ends before the steady state is attained* 

The low energy spectrum is very steep and process a cut-off (see 
Figures 2 and 3), which indicates the adiabatic deceleration is not very 
effective. The spectral index for accelerated particles (p > p 0 ) starts 
from values larger than -4 and then relaxes to pass -4 (see Figure 3)* 2 g 
In Figures 4 to 5 we used £<>= 3 x 10-^ erg and p =0*3, K|“ 10 

cm 2 sec\ ai^ao^O, and T «10 b K, nT*=0.1cm~, i.e. ^ = 0.572 and 1.79* 

The absorbing do undary is at a radius 5R. The spectrum (see Figures 4 
to 5) is very steep due to the large curvature of the shock and the 
escape of particles out of the system. 

We have studied the case of momentum dependent diffusion 
coefficient. The basic features are the same; some detailed changes 
are consistent with results of a steady state plane shock. 

V. Conclusion . The results discussed above may be understood as a 
consequence of the length scale or time scale involved (e.g. Prishchep 
and Ptuskin, 1980). There are three time scales involved: diffusion 

time scale Bi /K, convection time scale R/R and acceleration time scale 
K/R k . If the radius of curvature parameter RR/K is large, then the 
acceleration is very efficient and the spectrum will approach that from 
plane shock and vice versa. 
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